SMAD9-MYCN positive feedback loop represents a unique dependency for MYCN-amplified neuroblastoma

Background Neuroblastoma (NB) is the most common extracranial solid tumor occurring during childhood and high-risk NB patients have a poor prognosis. The amplified MYCN gene serves as an important determinant of a high risk of NB. Methods We performed an integrative screen using public NB tissue and cell line data, and identified that SMAD9 played an important role in high-risk NB. An investigation of the super-enhancers database (SEdb) and chromatin immunoprecipitation sequencing (ChIP-seq) dataset along with biological experiments of incorporating gene knockdown and CRISPR interference (CRISPRi) were performed to identify upstream regulatory mechanism of SMAD9. Gene knockdown and rescue, quantitative real-time PCR (Q-RT-PCR), cell titer Glo assays, colony formation assays, a subcutaneous xenograft model and immunohistochemistry were used to determine the functional role of SMAD9 in NB. An integrative analysis of ChIP-seq data with the validation of CRISPRi and dual-luciferase reporter assays and RNA sequencing (RNA-seq) data with Q-RT-PCR validation was conducted to analyze the downstream regulatory mechanism of SMAD9. Results High expression of SMAD9 was specifically induced by the transcription factors including MYCN, PHOX2B, GATA3 and HAND2 at the enhancer region. Genetic suppression of SMAD9 inhibited MYCN-amplified NB cell proliferation and tumorigenicity both in vitro and in vivo. Further studies revealed that SMAD9 bound to the MYCN promoter and transcriptionally regulate MYCN expression, with MYCN reciprocally binding to the SMAD9 enhancer and transactivating SMAD9, thus forming a positive feedback loop along with the MYCN-associated cancer cell cycle. Conclusion This study delineates that SMAD9 forms a positive transcriptional feedback loop with MYCN and represents a unique tumor-dependency for MYCN-amplified neuroblastoma. Supplementary Information The online version contains supplementary material available at 10.1186/s13046-022-02563-3.


Introduction
Neuroblastoma (NB) is the most common extracranial solid tumor occurring during childhood. NB threatens children's quality of life and survival, and accounts for 8-10% of childhood cancers and 15% of cancer-related deaths in children [1,2]. Based on age, the International Staging System (INSS), histology, genomic abnormalities, MYCN amplification and other characteristics, NB patients are classified into very low-, low-, intermediate-and high-risk groups according to the Children's Oncology Group (COG) [3,4]. DNA amplification of the MYCN proto-oncogene is an important determinant of a high risk in NB patients [3,4]. Despite multimodal treatment strategies, the overall survival (OS) rate for high-risk NB patients remains below 60% [1,5,6].
In addition to MYCN amplification, gene abnormalities including ALK mutation [7], MYC amplification [8], TERT rearrangement [9], and ATRX inactivation [10], are risk factors for NB. However, > 25% of high-risk NB patients do not show the abovementioned genetic abnormalities [9,11]. Moreover, NB is a specific tumor characterized by transcriptional gene abnormalities and MYCN is the most critical transcription factor (TF) [12,13]. Therefore, more transcriptional biomarkers must be identified to improve the clinical evaluation of high-risk NB patients.
Most recent clinical trials for NB are considered "basket trials", which are trials utilizing targeted drugs that are universally applicable for many cancer types [14,15]. Nevertheless, neuron-derived NB is a specific type of pediatric cancer characterized by epigenetic and transcriptional aberrations [16][17][18]. To illustrate this uniqueness, NB is a super enhancers (SEs)-driven tumor characterized by high expression of SEs-targeted genes such as MYCN, PHOX2B, GATA3, HAND2, TBX2, ISL1, ASCL1, and LMO1 in transcriptional core regulatory circuitries (CRCs) [13,17,19]. Here we performed an integrative analysis of public datasets and found that SMAD9 was characterized by specific high expression and dependence among all cancer types and was positively correlated with MYCN expression and a poor prognosis in high-risk NB patients. In our models, SMAD9 promoted MYCN-amplified NB cell growth in vitro and in vivo, and its upstream and downstream targets were strongly associated with MYCN. High SMAD9 expression was associated with MYCN-associated autonomous neural tumorigenicity and a robust cancer cell cycle.

Materials and methods
The study was approved by the Institutional Review Board (IRB) of Shanghai Children's Hospital (SCH), Shanghai Jiao Tong University, in accordance with the principles of Declaration of Helsinki. Written informed consent to participate in this study was provided by the participants' legal guardian or next of kin. Patients' identities and privacy were protected and anonymized in the study. In addition, the Medical Experimental Animal Administrative Committee in Shanghai approved all animal experimental protocols.

Patients and specimens
Patients were included based on the criteria for highrisk NB reported by the COG [3]: 1) MYCN amplification despite any other conditions; 2) INSS = 4 and age ≥ 1.5 years; and 3) INSS = 4 and the existence of genomic abnormalities. Finally, 14 high-risk NB samples obtained between January 2015 and December 2019 were selected for clinical validation. These tissue samples were excised in the operating room, immediately snap-frozen in liquid nitrogen and then stored at − 80 °C.

Quantitative reverse transcription PCR (Q-RT-PCR)
Approximately 25 mg tissue samples or 2.5 × 10 5 cells were lysed in TRIzol reagent (Thermo Fisher Scientific, #TR118) for total RNA extraction. RNA was reverse transcribed to cDNA with a High-Capacity RNA-to-cDNA kit (Thermo Fisher Scientific, #4387406). Quantitative PCR was performed using an Applied Biosystems Quant-Studio ™ 5 Real-Time PCR System (Thermo Fisher Scientific, #A34322) with SYBR Green Master buffer (ROX) (Thermo Fisher Scientific, #A25742). GAPDH was used as an internal control, and mRNA levels were calculated with the 2^delta Ct approach. The primer sequences are summarized in Supplemental Table 3.

Cell viability assay
Cells were plated in 96-well plates in triplicate at a density of 750 cells/well (small cells except for SK-N-SH and RPE-1 cells) or 375 cells/well (large cells) in 100 μL of culture medium. Then the CellTiter-Glo ® luminescent cell viability assay (Promega, #G7573) was performed to assess cell viability at days 0, 2, 4 and 6 according to the manufacturer's protocol.

Clone formation assay
A total of 375 cells/well (small cells) or 187 cells/well (large cells) were seeded in 12-well plates. Fresh medium was added every 5 days in the first 7 days. Half of the medium was removed, and the same amount of fresh medium was added every 4 days afterward. Generally, colonies were visible after 8-21 days, and the cells were washed with PBS, fixed with a 10% neutralized formaldehyde solution, and then stained with 0.5% crystal violet (Sigma, #C6158-100G) containing 25% methanol.

SMAD9 constructs
The sequence of the SMAD9 coding sequence (CDS) was screened to design a spot mutation without a change in amino acid sequence in the shSMAD9 target site. Wild-type (WT) and mutated (MUT) SMAD9 CDSs were cloned into the G418-resistant vector plenti-EF1α-MCS-IRES-Neo-WPRE (EcoRI and BamHI) by homologous recombination using the ClonExpress kit (Vazyme, #C112). The empty vector (EV), the rebuilt SMAD9-WT and SMAD9-MUT plasmids were applied for virus packaging in HEK293T cells and NB stably transfected cells (STCs) were established with G418 selection. Thereafter, STCs were transfected with shSMAD9 and selected with puromycin.

NB patient-derived Xenograft (PDX) and cells (PDCs)
Two primary NB samples from high-risk patients (SCH-NB010# and SCH-NB016#) were orthotopically xenografted into NSG mice in the first passage. Then, xenografts in the second, third and fourth passages in vivo were subcutaneously transplanted (passaged) into NSG mice. PDCs were extracted in the fourth passage and cultured in DMEM/F12 supplemented with 10% FBS. According to previous reports [22,23], the formation of neuron-shaped cells indicated good cell line generation. After the formation of large numbers of neuron-like cells, NB PDCs were digested and transfected with shSMAD9 virus for further analysis.

Tetracycline-inducible SMAD9 knockdown (Tet-on)
The shSMAD9 sequences were cloned into the Tet-on puromycin-resistant plasmid (Addgene, #21915). Empty vector (EV) and Tet-on-shSMAD9 plasmids were applied for virus packaging in HEK293T cells and for NB STC establishment with puromycin selection. Doxycycline (Dox) was used in Tet-on STCs. Therefore, STCs were expanded for validation in vitro (1 μg/mL Dox to induce SMAD9 knockdown) and tumor growth in vivo (2 mg/ mL Dox with 2% sucrose).

Animal experiments
BALB/c nude female mice (4-6 weeks) were purchased from the Experimental Animal Center of the Chinese Academy of Sciences (Shanghai, China). For subcutaneous cell line xenografts, 5 × 10 6 NB Tet-on STCs were subcutaneously transplanted in the dorsal flanks of mice on each side. When the tumors reached the volume of approximately 100 mm 3 , the mice were randomly divided into two groups and provided water containing 2 mg/mL Dox with 2% sucrose or 2% sucrose as a control. Tumor volumes were measured every 3 days with the formula 1/2 (long axis * short axis^2). Mice with tumors > 1500 mm 3 were euthanized.

Histological analysis
Histological analyses, including hematoxylin-eosin (HE) and immunohistochemical (IHC) staining, were performed by Servicebio Biotechnology Company (Shanghai, China). IHC staining was performed using primary antibodies against Ki-67 (Servicebio, #GB121499) and cleaved caspase-3 (Cell Signaling Technology, #9661). The stained cells were counted as the percentage of total positive cells in the five random fields of view using the IHC profiler plugin in the ImageJ software (v1.52p, USA).

ChIP-seq and ChIP-qPCR analysis
The Flag-tagged SMAD9 CDS was reconstructed in the plasmid pCDH-CMV-MCS-EF1-Puro using the Clon-Express kit (Vazyme, #C112). Empty vector (EV) and Flag-tagged SMAD9 CDS plasmids were applied for virus packaging in HEK293T cells and NB STC establishment with puromycin selection. Thereafter, NB STCs (1 × 10 7 ) were digested with micrococcal nuclease (New England Biolabs, #M0247S) and sonicated for ChIP experiments. Briefly, 10 μg of chromatin were immunoprecipitated with a Flag antibody (Cell Signaling Technology, #14793S), or H3K27Ac antibody (AM39133, Active Motif ) and 50 μL Pierce ChIP-grade Protein A/G Magnetic Beads (Thermo Scientific, #26162). Immunoprecipitated DNA and INPUT DNA were then purified and sequenced by Romics (Shanghai, China). ChIP-seq data were first mapped against the human genome build hg19 using BowTie2 [24]. Then, a model-based analysis of ChIP-seq data (MACS2: v2.2.7.1) was used to identify peak regions of ChIP-seq enrichment [25]. Peaks were annotated by performing a HOMER analysis [26]. Deep-Tools was used to generate BigWig files for visualization [27]. ChIP-seq data from biological experiments and public datasets were subjected to the same pipeline.
For the ChIP-qPCR analysis, 5-10 × 10 6 cells were harvested and 5-10 μg of chromatin were immunoprecipitated with a Flag antibody and 50 μL of Pierce ™ ChIP-grade Protein A/G Magnetic Beads. qPCR was performed using SYBR Green buffer (Thermo Fisher Scientific, #A25742). The primers were designed according to the peak sequence of Flag on the promoter of MYCN (Supplemental Table 4).

Luciferase reporter assay
A luciferase reporter assay was performed as described in a previous study [28]. The MYCN promoter or SMAD9 enhancer sequences were cloned into pGL3-luc vectors according to the peak sequences of SMAD9-Flag to MYCN or MYCN to SMAD9 based on ChIP-seq data.
The cloning PCR primers are presented in (Supplementary Table 6).
For MYCN knockdown, we constructed BE(2)-C Teton shMYCN cells by infecting wild-type BE(2)-C cells with the pLKO-Tet-On-shMYCN virus followed by continuous puromycin selection. Then, BE(2)-C Tet-on shMYCN cells were transfected with pGL3-luc based reporter plasmids and Renilla pGL3-Rluc control plasmid (ratio, 10:1). After 24 hours, transfection media were replaced with media containing or lacking doxycycline (1 μg/ml), and cells were harvested 48 h later. Procedures similar to those described above were used for BE(2)-C Tet-on shSMAD9 cells. Samples were assayed with a dual luciferase assay system (Promega, E1910). Firefly luciferase activity was normalized to Renilla luciferase activity.

Analysis of RNA sequencing data
For the RNA sequencing analysis of our own samples, trim galore was used to automatically detect and trim adaptors. Reads were mapped to the hg19 reference genome using HISAT2. Read counts were generated using HTSeq (version 0.11.1) and fragments per kilobase million (FPKM) values were calculated from the number of reads that mapped to each gene sequence, and the gene length and sequencing depth were considered.

Identification of differentially expressed genes (DEGs) and construction of a Venn diagram
The DEGs identified among NB tissues vs. normal adrenal glands, NB cell lines vs. neural crest cells and shScramble vs. shSMAD9 MYCN-amplified cells were analyzed using the false discovery rate (FDR) moderated limma test (package "DESeq2" in R). The cutoff was set for DEG selection based on the criterion |log 2 (fold change, FC)| > 1 with P < 0.05 in public datasets and |log 2 FC| > 0.4 with P < 0.05 in shSMAD9 RNA-seq data. A Venn diagram was constructed to visualize the overlapping genes using an available online tool (http:// bioin forma tics. psb. ugent. be/ webto ols/ Venn/).
In addition, we extracted the top 500 genes (ranked by P value) after MYCN knockdown (KD) from GSE80397 and GSE121529 and built the gene sets MYCN_TAR-GETS_GSE80397 and MYCN_TARGETS_GSE121529. An FDR q-value < 0.25 was considered statistically significant.

Protein-protein interaction (PPI) network, Cytoscape and Gene Ontology (GO)/Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses
The STRING (http:// string-db. org/) database was applied to determine the PPI network of 784 common downregulated genes with an interaction score of 0.4, after which Cytoscape combined with the CytoHubba plugin was used to visualize the PPI network and hub genes. GO/KEGG enrichment analyses were performed using the "clusterProfiler" and "org.Hs.eg.db" packages in R, and the top items were visualized.

Statistical analyses
GraphPad Prism v9.2.0 or R v3.6.3 software was applied for statistical analysis. Comparisons between two groups were performed using unpaired two-tailed Student's t test, and comparisons among more than two groups were performed using one-way ANOVA. Tumorigenicity in vivo was compared between the two groups using twoway ANOVA. The log-rank (Mantel Cox) test was used to analyze the Kaplan-Meier survival curves. P < 0.05 was considered statistically significant.

High SMAD9 expression is an indicator of a poor prognosis for a subset of high-risk patients with NB
The workflow of the whole study is shown in Fig. S1. An integrative analysis was performed by utilizing a DepMap dataset of the 1st ranked highly expressed and dependent genes (Figs. S2A and S2B) and public high-risk NB tissue and cell line datasets (Figs. 1A and S2A) to identify transcripts that were specifically expressed in high-risk patients with NB. We identified 16 potential NB riskrelated genes (Fig. 1A) and ranked them by fold change in expression (NB compared with other cancer types; Fig.  S2C). Notably, we found that SMAD9 expression levels were more than 6-fold higher in NB, similar to the wellknown transcripts including PHOX2A, MYCN, TBX2 and GATA3, indicating specificity (Figs. 1A and S2C). In NB, SMAD9 exhibited high expression among all cancer and normal cells (Figs. 1B and S2D). Moreover, SMAD9 expression was positively associated with high expression of MYCN and the MYCN amplification status in cell lines (Figs. 1B and S2D). An analysis of public functional genomics data showed that NB tumor viability depended on SMAD9, and its dependency score was potentially correlated with the MYCN amplification status and MYCN expression level (Fig. 1C). In addition, we detected a higher SMAD9 expression in high-risk NB tissues compared with normal adrenal gland (AG) tissues, and SMAD9 expression was higher in MYCN-amp and MYCN-high NB tumors (Fig. 1D). Interestingly, lowstage or nonhigh-risk NB tumors exhibited relatively higher expression of SMAD9 (Fig. S2E), therefore, the correlation between SMAD9 and MYCN expression was significantly diminished in nonhigh-risk NB (Fig. S2F) and all NB tumors (Fig. S2G).
We selected several public datasets and collected clinical information from high-risk NB patients in Shanghai Children's Hospital (SCH) to further analyze the prognostic significance of SMAD9 (Fig. S3A). High SMAD9 expression indicated a poor prognosis for patients in the high-risk cohorts (Figs. 1E and S3B-S3C). Similar to the results of SMAD9 expression in NB tumors (Figs. 1D and S2E-S2G), SMAD9 was not an indicator of inferior outcome in nonhigh-risk NB (Fig. S3D) and all NB patients (Fig. S3E). MYCN-amp NB was a subset of high-risk NB, and high SMAD9 expression also indicated a poor prognosis for MYCN-amp NB patients (Figs. S3C and S3F). Notably, high-risk NB tissue samples with high expression of SMAD9 were enriched in gene sets of GSEA concerning oncogenic malignancy terms (Fig. 1F). Taken together, our integrative analyses indicate that SMAD9 is expressed at high levels in NB cell lines and tissues and is associated with a poor prognosis of high-risk patients.

SEs-targeted SMAD9 is potentially activated by NB-specific genes in CRC
We performed an upstream regulatory analysis in SEdb to investigate the reason why SMAD9 was highly and specifically expressed in NB, showing that SMAD9 was a specific SEs-targeted gene in NB compared to other cancer types according to SE counts ( Figs. 2A and S4A-B), SE rank (Fig. 2B) and SE binding patterns (Fig. 2C). In addition, NB is a peripheral nervous tumor with high specificity driven by SE-targeted genes in CRC [13,17]. Therefore, we analyzed public ChIP-seq data (GSE94822) and found that MYCN, GATA3, PHOX2B and HAND2 (TFs in CRC) shared peaks with H3K27ac peaks at the enhancer region of SMAD9 in MYCN-amplified BE(2)-C cells (Fig. 2D). Data from the DepMap portal showed that SMAD9 expression was positively correlated with these SE-targeted genes (Fig. 2E).
We investigated some knockdown datasets to validate the binding potency of these SEs-targeted genes on SMAD9 enhancer region and found that SMAD9 expression was downregulated by knocking down MYCN (GSE121529) and PHOX2B (GSE124451) (Fig. S4C). In addition, we used shRNAs to knockdown MYCN, GATA3, PHOX2B and HAND2. Interference with one of these 4 genes decreased SMAD9 mRNA levels in BE(2)-C cells (Fig. 2F) and SK-N-BE2 cells (Fig. 2G). Moreover, the CRISPRi assay showed that the interruption of the binding of CRC signals to SMAD9 reduced SMAD9 expression in NB cells (Fig. 2H). The disturbance of the binding sites of the CRC signal to SMAD9 reduced NB cell viability (Fig. 2I), and colony formation (Fig. 2J). Taken together, these data suggest that NB-specific TFs in CRC might bind the SMAD9 super-enhancers region and induce high SMAD9 expression in NB.

SMAD9 knockdown suppresses MYCN-amplified NB growth
We detected MYCN expression in NB cell lines (SK-N-SH, SH-SY5Y, SK-N-BE2, BE(2)-C, IMR-32) at the protein level (Fig. S5A) and SMAD9 and MYCN mRNA expression in NB cell lines, a fibroblast cell line (RPE-1), epithelium-derived cell line (MCF-10) and NB tissues from SCH (Fig. S5B). We did not observe a significant difference in SMAD9 expression between MYCN-nonamplified and MYCN-amplified cell lines and tissues, but MYCN-high NB tissues (classified by the median expression) had a significantly higher SMAD9 level than MYCN-low NB tissues (Fig. S5B). In addition, we analyzed SMAD9 expression in noncancerous cells using Q-RT-PCR (Fig. S5B) and in datasets ( Fig. S2D and S5C), and observed lower expression of SMAD9 in noncancerous cells. Next, we used shRNAs to knock down SMAD9 expression in MYCN-amplified cells (BE2C, IMR32 and SK-N-BE2 cells; Fig. 3A), and SMAD9 knockdown significantly inhibited NB growth in vitro (Fig. 3B). In addition, a colony formation assay indicated that SMAD9 inhibition suppressed MYCN-amplified NB growth (Fig. 3C and D). However, SMAD9 knockdown in MYCN nonamplified cells (SK-N-SH and SH-SY5Y cells; Fig. 3E) did not significantly repress NB growth in vitro (Fig. 3F) or even promote NB growth (lower panel; SH-SY5Y cells in Fig. 3F). Similarly, SMAD9 knockdown in MYCN nonamplified cells did not inhibit colony formation (Fig. 3G and H). These results were consistent with functional genomics data in the DepMap portal (Fig. S5D). Moreover, we generated two types of MYCN-amplified patient derived cells (PDCs named SCH-NB010# and SCH-NB016#) and found that knockdown of SMAD9 suppressed NB growth (Figs. S5E and S5F). Finally, reconstitution of SMAD9 expression (Fig. 3I) increased the viability of MYCN-amplified NB cell lines (Fig. 3J). Based on these results, SMAD9 knockdown inhibits MYCN-amplified NB growth in vitro.
We employed Tet-on shRNA methods to induce SMAD9 knockdown in BE(2)-C Tet-on-shSMAD9 STCs as a method to validate the function of SMAD9 in MYCN-amplified NB in vivo (Fig. S6A), and found that doxycycline-induced SMAD9 knockdown suppressed NB growth in vitro (Figs. S6B-S6D). Afterward, we suppressed SMAD9 expression in vivo and found that SMAD9 knockdown in NB xenografts led to a decreased tumor volume and mass ( Fig. 4A-C). SMAD9 knockdown was validated using a Q-RT-PCR assay (Fig. 4D). The histological analysis indicated that SMAD9 knockdown abolished the tumor microstructure (Figs. 4E and S6E), decreased cell proliferation and increased apoptosis (Figs. 4E, F and S6E). Collectively, the suppression of SMAD9 attenuates MYCN-amplified NB growth both in vitro and in vivo.

The SMAD9 binding pattern is associated with MYCN
Next, we investigated the transcriptional role of SMAD9 (Fig. S7A). Due to the lack of commercially reliable antihuman SMAD9 antibodies, we constructed SMAD9 overexpression constructs with the 3xFlag tag plasmid and established stably transfected MYCN-amplified BE(2)-C cells (Figs. S7B and S7C). We conducted ChIP-seq for BE(2)-C-SMAD9-Flag STCs, combined the results of RNA-seq for BE(2)-C-shSMAD9 cells and found 7 overlapping genes (Fig. 5A). Q-RT-PCR analyses validated that these genes were potential downstream genes (Fig. 5B). Based on the expression and dependency data from the DepMap portal, we found that MYCN was the only gene that exhibited a high and specific dependency on NB (Fig. 5C).
Specifically, we visualized strong binding for SMAD9 and H3K27ac in the MYCN promoter region (Fig. 5D). ChIP-QPCR experiments validated the aforementioned finding (Fig. 5E). CRISPRi experiments showed that blocking SMAD9 binding sites suppressed MYCN mRNA (Fig. 5F) and protein expression (Fig. 5G) in three types of MYCN-amplified dCas9 STCs. In parallel, the disruption of SMAD9 binding sites for MYCN diminished cell proliferation (Fig. 5H) and inhibited colony formation in MYCN-amplified dCas9 STCs (Figs. 5I and S7D). Furthermore, the recovery of SMAD9 expression restored MYCN mRNA (Fig. 5J) and protein expression (Fig. 5K). We further validated the transcriptional feedback loop between SMAD9 and MYCN by constructing a luciferase reporter containing the MYCN promoter or SMAD9 enhancer. The dual-luciferase reporter analysis showed that SMAD9 or MYCN knockdown reduced the activity of MYCN promoter or SMAD9 enhancer (Fig. 5L).
Interestingly, we also detected a strong binding signal for SMAD9 and H3K27ac in the MYCNOS promoter region (Fig. S7D). MYCNOS expression was higher in high-risk NB tissues than in nonhigh-risk NB tissues, and higher in MYCN-amp NB tissues than in nonMYCNamp NB tissues (Figs. S7E and S7F). Disrupted binding sites of SMAD9 to MYCN also suppressed MYCNOS transcription (Fig. S7G) and the recovery of SMAD9 expression restored MYCNOS mRNA levels (Fig. S7H). Similar to the results from ChIP-seq, the GESA of RNAseq data from BE(2)-C and IMR-32 cells showed that high SMAD9 expression had a high normalized enrichment score in the terms of "WEI_MYCN_TARGETS_ WITH_E_BOX", "KIM_MYCN_AMPLIFICATION_UP", "MYCN_TARGETS_GSE121529", "MYCN_TARGETS_ GSE80397" and other MYCN-associated signatures ( Fig. 6A and B). In addition to NB cells, we observed a positive correlation between SMAD9 and MYCN (or MYCNOS) expression in high-risk NB tumors and a much lower correlation in the NB cohorts that included fewer high-risk patients (Fig. 6C), consistent with previous results (Figs. 1D and S2E-G). GSEA of NB tissues showed that high expression of SMAD9 had a high enrichment score in these MYCN-associated signatures (Fig. 6D).
Taken together, our results indicate that SMAD9 acts as a TF and is capable of binding MYCN to promote MYCN-amplified NB cell growth and that SMAD9-MYCN forms a positive transcriptional feedback loop.

Malignant cell cycle changes occur in response to SMAD9 knockdown
Finally, we performed unbiased transcriptomic analyses of shSCR and shSMAD9 samples among BE(2)-C and IMR-32 cells, with > 25% decreases in MKI67, BCL2 and CCND1 levels (Fig. 7A). GSEA showed that high SMAD9 expression was mainly enriched in terms related to the malignant cell cycle both in BE(2)-C and IMR-32 cells (Figs. 7B and S8A-B). Therefore, we combined the RNA-seq results from the two cell lines and subsequently investigated the 784 overlapping downregulated genes (Fig. 7C). GO/KEGG analyses revealed that the 784 genes were mainly enriched in cancer cell cycle terms such as "chromosome segregation" and "mitotic nuclear division" in GO-biological processes (GO-BP) analysis (Fig. 7D). Similarly, other signatures, such as "condensed chromosome" in GO cellular component (GO-CC) terms, "microtubule binding" in GO molecular function (GO-MF) terms and "cell cycle" in KEGG pathways, revealed the significance of phenotypic malignancy (Fig. S8C). Regarding upregulated genes in shSMAD9 cells, we identified 570 common upregulated genes, and GO/KEGG analyses suggested the enrichment of mature nervous signatures such as "modulation of chemical synaptic transmission" in GO-BP terms and "neuron to neuron synapse" in GO-CC terms (Fig. S8D).
The common downregulated genes were subjected to PPI analysis to further identify the downregulated hub genes, and we found that AURKA-linked cell cycle and mitotic nuclear division genes such as AURKA, AURKB and CCNA2 might play a pivotal role in SMAD9-MYCNmediated malignancy (Fig. 7E), which was validated by performing Q-RT-PCR analysis of NB cells (Fig. 7F) and tissues from SCH (Fig. 7G). We briefly validated previous findings concerning the MYCN and cell cycle pathways [30][31][32] by analyzing two datasets from MYCN knockdown cells and found that most of SMAD9-high hub genes had comparable edges in MYCN-high NB cells (Figs. S9A and S9B). The validated cell cycle genes were downregulated in MYCN-low NB cells and tumors (Fig.  S9C). Notably, the ChIP-seq data for MYCN suggested the significance of E-box motifs (CAC GTG or CANNTG) in MYCN-amp NB cells (Figs. S9D and S9E). We identified abundant E-box sequences in the promoters of the hub cell cycle genes in different types of MYCN-amp NB cells (Figs. S9F-S9G).
The combined analysis of SMAD9 binding patterns and transcriptomic changes in SMAD9 knockdown cells indicated that SMAD9 might contribute to MYCN-mediated autonomous tumorigenicity and the cancer cell cycle in MYCN-amplified NB. As described previously, SMAD9 is a NB-specific SEs-targeted gene mediated by MYCN and other TFs in CRC and dual luciferase assays showed that SMAD9 or MYCN knockdown reduced the activity of MYCN promoter or SMAD9 enhancer. Collectively, the positive feedback loop between SMAD9 and MYCN reinforces autonomous nervous system tumorigenicity and the tumor cancer cycle in MYCN-amplified NBs (Fig. 7H).

Discussion
Despite the use of multimodal therapeutic strategies, COG-defined high-risk NB has a low 5-year OS rate [33]. NB is driven by the interplay of regulatory TFs and dynamic alterations in chromatin structure [34]. Therefore, in this study, we performed integrative analyses of public datasets and identified SMAD9 as a potential prognostic risk factor in the high-risk NB patients. In addition, SMAD9 was shown to be a SEs-targeted gene in NB, and we revealed SMAD9 as a regulator of MYCNamplified NB cell proliferation and tumorigenicity. Mechanistic studies showed that SMAD9 bound to the MYCN promoter, partially regulated MYCN expression and further mediated the cancer cell cycle. Taken together, MYCN, GATA3, PHOX2B and HAND2 reciprocally triggered SMAD9 transcription by binding its enhancers, thus forming a positive feedback loop between SMAD9 and MYCN.
MYCN has been implicated in NB development [30], drug resistance [45] and metastasis [46]. Recent studies have aimed to identify targets that indirectly suppress MYCN expression [21,34,[47][48][49]. To our knowledge, this study is the first to identify a direct link between SMAD9 and MYCN. First, we found that SMAD9 had a specific expression and dependency profile in NB and showed a positive correlation with MYCN expression. Next, we showed that MYCN coordinated with other TFs in CRC to target SMAD9 enhancers and transcriptionally activate SMAD9 expression. In addition, SEs-targeted SMAD9 promoted the growth of MYCN-amplified NB cells both in vitro and in vivo. Furthermore, a combination of ChIP-seq, RNA-seq, and CRISPRi analyses indicated that SMAD9 targeted the MYCN promoter and transcriptionally activated MYCN expression. Along with the dual-luciferase reporter assays, our results highlight the significance of the positive feedback loop between SMAD9 and MYCN. Some of the existing positive feedback loops between other tumor-dependent genes and MYCN have been reported [21,34], and our study has provided additional evidence for MYCN networks.
It has been documented that MYCNOS functions as an antisense RNA and results in a reduction of upstream MYCN promoter usage and thus upregulates MYCN expression [50]. In addition, CCCTC-binding factor (CTCF) cooperates with noncoding RNA MYCNOS to promote neuroblastoma progression through facilitating MYCN expression [51]. We identified a strong binding for SMAD9 and H3K27ac in the MYCNOS promoter region and the recovery of SMAD9 increased MYCNOS expression, suggesting that MYCNOS might be critical in the SMAD9-MYCN positive feedback loop.
MYCN induces quiescent cells to reenter the cell cycle and shortens the cell cycle time, in particular shortening G1 phase [30,31,[52][53][54]. SMAD9 knockdown attenuated the transcriptomic phenotypes of MYCNassociated autonomous nervous system development and the cancer cell cycle. We identified overlapping cell cycle hub genes, such as AURKA and CCNA2, in both MYCN and SMAD9 knockdown NB cells. And we validated the E-box sequences in the promoter of the hub cell cycle genes. Taken together, these results suggest that SMAD9 potentially mediates the MYCN-relevant cell cycle in a subset of high-risk NB cells.
Some limitations of our study should be noted. We only used BE(2)-C cells for in vivo experiments because other NB cells exhibited less tumorigenic potential and high heterogeneity in nude mice. Although we identified a possible positive feedback loop between SMAD9 and MYCN, our current study did not precisely delineate the crosstalk between the SMAD9 pathway and MYCN expression. Notably, due to the lack of a reliable antihuman SMAD9 antibody, we constructed a SMAD9-Flag-tagged vector in MYCN-amplified cells to perform ChIP-seq analysis according to previous reports [55,56]. Although our results revealed a low number of SMAD9-Flag binding peaks, we did not obtain evidence for false positive or negative binding events in Flag ChIP-seq experiments [57].
In summary, high SMAD9 expression in NB was specifically regulated by multiple TFs at the SEs region in CRC. In addition, SMAD9 knockdown exerted an antitumor effect on MYCN-amplified NB by disrupting the positive feedback loop between SMAD9 and MYCN. These findings will help develop novel diagnostic and therapeutic strategies for high-risk NB.

Conclusion
SMAD9 forms a positive transcriptional feedback loop with MYCN and represents a unique tumor dependency for MYCN-amplified neuroblastoma.